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Abstract 


We present a new and accurate method to render the atmosphere in real time from any viewpoint from ground 
level to outer space, while taking Rayleigh and Mie multiple scattering into account. Our method reproduces 
many effects of the scattering of light, such as the daylight and twilight sky color and aerial perspective for all 
view and light directions, or the Earth and mountain shadows (light shafts) inside the atmosphere. Our method is 
based on a formulation of the light transport equation that is precomputable for all view points, view directions and 
sun directions. We show how to store this data compactly and propose a GPU compliant algorithm to precompute 
it in a few seconds. This precomputed data allows us to evaluate at runtime the light transport equation in constant 
time, without any sampling, while taking into account the ground for shadows and light shafts. 


Categories and Subject Descriptors (according to ACM CCS): 1.3.7 [Computer Graphics]: Three-Dimensional 


Graphics and Realism 


1. Introduction 


Atmospheric effects are very important to increase the real- 
ism of outdoor scenes in many applications. The sky color 
gives key indications about the hour of the day, and the 
aerial perspective gives an important cue to evaluate dis- 
tances. Rendering these effects in real time, continuously 
from ground to space, is desirable in many games or applica- 
tions, such as flight simulators or Earth browsers like Google 
Earth. This is especially true for applications that target re- 
alism, such as Celestia or Nasa WorldWind. However, these 
applications currently use very basic models to render these 
effects, which do not give realistic images. 


In this paper we propose a method to render these effects 
in real time, from any viewpoint from ground to space. This 
method accounts for multiple scattering, which is important 
to correctly render twilight, or the shadow of the Earth in- 
side the atmosphere (see Figure 8). It is based on moderate 
simplifying assumptions that allow us to get a better approx- 
imate solution of the rendering equation (compared to pre- 
vious work), in which most terms can be precomputed. Our 
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method is the first real-time method accounting for all view- 
points, all view and sun directions, and multiple scattering. 


The next sections are organized as follows. Section 2 in- 
troduces the physical model and the rendering equation and 
reviews the related work. Section 3 presents our resolution 
method to get a precomputable formulation. Sections 4 and 5 
present our precomputation and rendering algorithms. Sec- 
tion 6 gives implementation details and presents our results. 


2. Atmospheric models 


Rendering atmosphere illumination relies on two aspects: a 
physical model of the local medium properties, and a simu- 
lation of the global illumination exchanges up to the viewer 
eyes. This includes exchanges with the ground, which can 
be modeled as a Lambertian surface with a height field of 
reflectance (x, À), normal n(x), etc. 


Most computer graphics (CG) papers, starting with 
[NSTN93], rely on a physical model of the medium com- 
prising air molecules and aerosol particles, summarized in 
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Section 2.1. However, the classical rendering equation for 
participating media is rarely completely accounted for in at- 
mospheric CG models, especially for interactive rendering. 
We restate the general model in Section 2.2 and we present 
its approximations in previous CG models in Section 2.3. 


2.1. Physical model 


The physical model commonly used in CG is a clear sky 
model based on two constituents, air molecules and aerosol 
particles, in a thin spherical layer of decreasing density be- 
tween Rg = 6360 km and R; = 6420 km (see Figure 1). 


At each point, the proportion of light that is scattered 6 de- 
grees away from its incident direction is given by the prod- 
uct of a scattering coefficient B° and of a phase function P. 
B° depends on the particule density and P describes the an- 
gular dependency. For air molecules B* and P are given by 
the Rayleigh theory [TS99]: 


Br(h,d) = 


Pals) = a= +R) 


where h = r — Rg is the altitude, À the wavelength, n the 
index of refraction of air, N the molecular density at sea 
level Rg, and Hr = 8 km is the thickness of the atmo- 
sphere if its density were uniform. As in [REK*04], we use 
By = (5.8, 13.5,33.1)10~° m7! for A = (680, 550,440) nm. 
Aerosols also have an exponentially decreasing density, with 
a smaller height scale Hy © 1.2 km. Their phase function is 
given by the Mie theory, approximated with the Cornette- 
Shanks phase function [TS99]: 


FR (D) 


where u = cos 0 (2) 


Bia (h, À) = Bir(0,A)e Fz (3) 


3 (= @)(14 2) 
8m (2+ g7)(1 +8? —2gu)3/? 


Unlike air molecules, aerosols absorb a fraction of the inci- 
dent light. It is measured with an absorption coefficient B,. 
which gives the extinction coefficient Bi, = Biy + Bi, (see 
Figure 6 for typical values — BR = Bp for air molecules). 
Note that the variation of the index of refraction with al- 
titude causes a small bending of rays (less than 2 degrees 
[HMS05]). We ignore it for simplicity. 


Pu (u) = (4) 


2.2. Rendering equation 


We recall here the rendering equation in a participating 
medium, applied to the atmosphere. We note L(x, v,s) the 
radiance of light reaching x from direction v when the sun is 
in direction s, and x, (x, v) the extremity of the ray x +rv (see 
Figure 1). Note that x, is either on the ground or on the top 
atmosphere boundary r = R;. The transmittance T between 
Xo and x, the radiance Z of light reflected at xo, and the radi- 
ance J of light scattered at y in direction —v are defined as 


Figure 1: Our method. Left: the reference solution includes single- 
scattering (a) and multiple-scattering (b) integrated from x to Xo, all 
accounting for occlusion. Right: our approximation. Integration is 
done from x to Xs, ignoring occlusion (implicit via the use of Xs). 
(a) is unchanged. (b) is affected by ignoring occlusion of secondary 
scatters (this yields both positive and negative bias, and effect is 
small anyway). 


Figure 2: Definitions. (a) the atmospheric transparency T results 
from absorption and out scattered light. (b) T[L] is the light L re- 
flected at Xo. It is null on the top atmosphere boundary. (c) J [L] 
is the light L scattered at y in direction —v. (d) SL] is the light 
scattered towards x between x, and x, from any direction. 


follows (see Figure 2): 


T(x,Xo) = exp i 


L B sw) (5) 


iE{RM} 


T[L](xo,S) = Axe) [Lt @s)@n(x)de. o0 (6) 
TGV) = f E ORUD sdo O 
T ic{R,M 


Note that Z is null on the top atmosphere boundary. With 
these notations the rendering equation is [TS99]: 

L(x, v,s) = (Lo + RIL] +S[L])(x,v,s) (8) 

Lo(x,V,8) = T (X, Xo)Lsun, or 0 (9) 

RIL] (x, v,s) = T (x, xo )T[L] (xo,s) (10) 


Siv = f TaI aD 


where Lọ is the direct sunlight Lsun attenuated before reach- 
ing x by T'(x,x0). Lo is null if v Æ s, or if the sun is occluded 
by the terrain, i.e., if x, is on the ground. R[L] is the light re- 
flected at x, and also attenuated before reaching x, and S[L] 
is the inscattered light, i.e., the light scattered towards x be- 
tween x and xo (see Figure 2). 


2.3. Previous rendering methods 


Equation 8 is very complex to solve. Hence, many simpli- 
fying assumptions have been made in CG to find approxi- 
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mate solutions that are easier to compute (see [S1002] for a 
survey). Most real-time methods ignore multiple scattering. 
In this case Equation 8 reduces to L = Lo + R[Lo] + S[Lo]. 
However, even S[Lo]| is quite complex to solve. Some au- 
thors propose analytical solutions at the price of idealisa- 
tions: flat Earth with constant atmosphere density [HP02], 
or without Mie scattering [REK*04]. The flat Earth hypoth- 
esis limits them to observers on the ground. Otherwise, S [Lo] 
is generally computed by numerical integration [NSTN93], 
which can be done in real time using low sampling [O° N05]. 
A notable exception is [SFE07] who rely on precomputa- 
tions of this integral. However, in order to reduce the num- 
ber of parameters, they only take into account the view and 
sun zenith angles, and neglect the angle between the view 
and sun directions. Hence they cannot reproduce, e.g., the 
Earth’s shadow inside the atmosphere. 


Ignoring multiple scattering as above is acceptable for 
daylight but not for twilight [HMSO5]. This is because sun- 
light traverses much less atmosphere during the day than 
during sunset or sunrise. Hence, some authors propose meth- 
ods to account for multiple scattering. [PSS99] fit the results 
of a double scattering Monte-Carlo simulation with an ana- 
lytical model, but their model is only valid for an observer on 
the ground. [NDKY96] and [HMSOS] use volume radiosity 
algorithms to compute multiple scattering, but their methods 
are far from real-time (minutes to hours per image). 


In this paper we propose a new method to render the sky 
and the aerial perspective in real time, from all viewpoints 
from ground to space, while taking multiple scattering into 
account. It is inspired by [SFE07] and extends it with multi- 
ple scattering, with the previously ignored view-sun angle 
parameter, with a better parameterization for the precom- 
puted tables, and with a new method for light shafts. 


3. Our method 


For efficiency and realism, our goal is to precompute L as 
much as possible, with only minimal approximations. Our 
solution is based on an exact computation for zero and sin- 
gle scattering, and uses an approximation of occlusion ef- 
fects to compute multiple scattering. In fact we take the de- 
tailed ground shape into account for zero and single scatter- 
ing, in order to get correct ground colors, shadows and light 
shafts. But we approximate it with a perfect sphere of con- 
stant reflectance to compute multiple scattering, to allow for 
precomputation. 


Notations Before presenting our method we need some no- 
tations and auxiliary functions. We note L = Lo + (R +S) [L] 
the solution of Equation 8 for the case of a perfectly spher- 
ical ground of constant reflectance &. Lo, R, S, žo, I and 
so on are defined as before, but for this spherical ground. 
Note that thanks to the ground’s spherical symmetry x and v 
can be reduced to an altitude and a view zenith angle. Hence 
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functions of x, v,s such as L or S[L] can be reduced to func- 
tions of 4 parameters (2 for x, v and 2 for s). Note also that L 
(resp. L) can be expressed with a series in the linear opera- 
tors R and S (resp. R and S), where the i term corresponds 
to light reflected and/or scattered exactly i times: 


L=1o+(R+8)|[Lo] +(R+S)[(R+S)[Lo]] +... 
= Lo + Li + Lo +... = Lo + Lx (12) 


Zero and single scattering We compute Lo and R[Lo] ex- 
actly, during rendering. For this we use a shadowing algo- 
rithm to compute the sun occlusion (see Equation 9), and a 
precomputed table for the transmittance T, which depends 
on only 2 parameters (see Section 4). S [Lo] is more compli- 
cated. It is an integral between x and x, but, due to the occlu- 
sion term in Lo, the integrand is null at all points y that are in 
shadow (this is what gives light shafts). We suppose here that 
these points are between xs and Xo (see Figure | — the gen- 
eral case is discussed in Section 5). Then the integral can be 
reduced to the lit segment [x, xs]. Moreover, occlusion can be 
ignored since it is already accounted for via xs, i.e., Lo can 
be replaced with Lo. This shows that S[Lo] = f% T J [Lo]. 
By rewriting this as [¥ T [Lo] — i TJ [Lo], extending an 
idea introduced in [O’NO5] and reused in [SFEO7], we fi- 
nally get a formulation using precomputable functions of 2 
and 4 parameters, T and S[Lo]: 


S[Lol(x, v,s) = S[Lo](x, v, 8) — T(x, xs)S[Lo](xs,v,s) (13) 


Multiple scattering As shown above Lo and L4 can be com- 
puted exactly despite the occlusion. Unfortunately account- 
ing for occlusion in the other terms L2 +. . . = R[L:]+S[L;] 
is much more difficult. Hopefully, in this case the occlusion 
can be approximated. Indeed, multiple scattering effects are 
small compared to single scattering during the day, while 
the ground contribution is small when it is not directly lit 
by the sun. So we approximate occlusion effects in S[L:] 
by integrating the contribution of multiple scattering, com- 
puted without occlusion, between x and xs. This yields both 
positive and negative bias (see Figure 1). Mathematically, 
this approximation gives S[L.] ~ J% TJ [Ls]. We also ap- 
proximate occlusion effects in R[L+] with the ambient oc- 
clusion of an horizontal hemisphere due to the ground’s tan- 
gent plane, +. This gives R[L:] ~ R[L:] with: 


RIL] = T (x,Xo) 


oso) Pate ne) se à td) 


2 
E[L«](Xo,8) = LL 0)œ()de, or 0 (15) 


By using the same rewriting rule as for Equation 13, and by 
noting S[L] Ix = S[L](x, v,s), we finally get: 


L~Io+R[Lo] + RIL:] + S[E]\x—T(%,%s)S[L]}x, (16) 


where the first three terms can be quickly computed with the 
help of precomputed 2D tables for T and €[L.], and where 


S[L] can be precomputed in a 4D table. We now show how 
to precompute them, in tables of a reasonnable size. 
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Algorithm 4.1: PRECOMPUTE(A,yders) 


T(x,v) To, v) 

AE(x,s) — E[Lol(x,8) 

AS(x,v,s) — S[Lo ix, v,s) 

E(x,s)— 0 

S(x, v,s) — AS(x,v,s) 

for i — 1 to i < Norders 
AJ(x,v,8) — JIT AE + AS] (x,v,s) 
AE(x,s) — E[T&AE + AS] (x, s) = E[AS](x,s) 

do AS(x,v,s) S Je T(x, y)AJ (y, v,s)dy 

E(x, s) — E(x, s) +AE(x,s) 
S(x, v,s) — S(x, v,s) + AS(x, v,s) 


4. Precomputations 


We precompute T(x,%o(x,v)) for all x,v in a 2D table 
T(x,v). Due to spherical symmetry, T depends only on 
r= ||x|| and u = v.x/r [O’N05]. As [O’N05], we then use the 
identity T (x,y) = T(x, v)/T(y,v), with v = (y—x)/|ly—x]. 


We precompute €[L.] and $[L] in two tables E and S with 
an algorithm that computes each scattering order L; one after 
the other. This algorithm uses three intermediate tables AE, 
AS and AJ containing after each iteration i £[L;], S[L;] and 
J [Li]. AE and AS are added to the result tables E and S at 
the end of each iteration (R is computed with the identity 
RIL](x, v,s) = T(x, Xo) 2E[L](Ko,8) — see Algorithm 4.1). 


Angular precision Since S is a 4D table its size increases 
very quickly with resolution. So we can only use a lim- 
ited angular resolution for v. This poses a precision prob- 
lem, which is however limited to the strong forward Mie 
scattering. In order to solve it we separate the single Mie 
scattering term from all the others in S, so as to apply 
the phase function at runtime. For this we rewrite S[Z] as 

PuSm|Lo] + PR$R[Lo] + S[L;]. We then store Cy = Sy[Lo| 
and Cx = SR[Lo] + S[L«]/Pr separately, which requires 6 
values per entry in S. If necessary, for efficiency, this can be 
reduced to 4 values per entry by storing only the red com- 
ponent Cy, of Cm. In this case the other components can 
be approximated with a proportionality rule between Sm [Žo] 


and Sa[Lo], which gives Cy ~ Cx ct p at. 


Parameterization In order to store S[L] into S we need a 
mapping from (x,v,s) into table indices in [0, 1]*. A simple 
solution is to use r = ||x|| and the cosinus of the view zenith, 
sun zenith, and view sun angles, u = v.x/r, us = s.x/r and 
v = v.s (mapped linearly from [R¢, R;] x [—1,1}° to [0, 1]*. 


The problem of this parameterization is that it requires 
a very high resolution in u to get a good sampling for the 
aerial perspective. Consider for instance an observer near 
the ground looking horizontally, with a mountain at distance 
d (see Figure 3). The aerial perspective is given by Equa- 
tion 16 as S(x,v,s) — T(x,xs)S(xs,v,s). Then u = 0 for x, 


Figure 3: Viewing angle parameter. Left: using u gives artifacts. 
Right: using uu = do/dy or do/dy solves the problem (using 128 
values for u or uy in the precomputed sky radiance table S). 


6360 6370 6380 6390 6400 6410 6420 
1 T T T 


0.8 


0.6 


Figure 4: Parameterization. u,, Uj, Uu, as functions of r, u, lis. 


and d/Vr2 +d? for xs, which gives Au = 0.016 < 1 for 

= 100 km. This too small value gives visible artifacts (see 
Figure 3). In order to solve this problem we rely on a bet- 
ter parameterization. We replace u with uu, defined as the 
ratio between the distance do = ||Xo — x|| and the distance 
dh (resp. dy) from x to the horizon (resp. to the atmosphere 
boundary “behind” the horizon — see Figure 3). In the previ- 
ous example dy ~ (R? = R?) 1/2 for x and Xs, while do ~ dy 
for x and dy —d for xs, which gives Au, = 0.11 >> 0.016 
for d = 100 km. With this mapping 128 samples for uy are 
sufficient to avoid the above artifacts. 


Another problem is that S is discontinuous at the hori- 
zon, due to the discontinuity of the length of the viewing 
ray here. Hence a continuous mapping yields linear interpo- 
lations across this discontinuity, which causes artifacts. We 
solve this problem by ensuring that uy is itself discontinuous 
at the horizon (see Figure 4). Finally, we use an ad hoc non 
linear mapping for r and us, chosen so as to get a better pre- 
cision near the ground and for sun zenith angles near 90°. 
So our mapping from (x, v,s) into [0, 1]* is finally defined as 
follows: 


ur =p/H 

uu = 1/2 + (ru + VA)/(2p) 
1/2— (ru — VA+H?)/(2p + 2H) 

— 45-06) 1(1 — 9-3) 

uy =(1+v)/2 


if ru < 0 and A > 0 


otherwise 


with p = (r° —R2)"/?, H = (R?—R2)1/2, and A = r — p?. 
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ae 


Figure 5: Evaluation of l. Left: due to the false boundaries b and 
c, the computed length Az — An.zg = Zg — Za + Zc — Zp is larger than l. 
Clamping this value to zg — Zmin fixes the problem. Right: viewpoint 
in shadow. Using only the extruded edges, Xo would be seen as lit 


and | would be equal to 0 instead of zg — Znear. Projecting the back 
faces (dashed line) on the near plane solves the issue [HHLHO5]. 


5. Rendering 


In order to render the sky and the aerial perspective we evalu- 
ate Equation 16 at each pixel. Lo can be efficiently computed 
using T. Computing R[Lo] involves T, a(xo) and n(xo), and 
a shadow test to determine whether Xo is lit. Finally E and 
S are used to compute R[L:] and S[Z]. As in [SFEO7], x is 
the camera position or, if in space, the nearest intersection 
of the viewing ray with the atmosphere boundary. The only 
remaining non-trivial parameter is xs, which depends on the 
terrain shadows and gives light shafts. 


Most light shaft algorithms use sampling or slicing to per- 
form a numerical integration along the viewing ray, with a 
shadow map to find which samples are lit. Up to 100 sam- 
ples per ray must be used to eliminate the artifacts due to the 
discrete sampling [IJTNO7]. We propose here a new method 
inspired from shadow volumes [HHLHOS5]. It does not rely 
on numerical integration, and therefore does not suffer from 
these artifacts. We first show that an exact computation is 
possible but not adapted to the GPU. We then present an ap- 
proximate solution better adapted to the GPU. Our idea is 
to use the precomputed integral S to compute the inscattered 
light due to each lit segment [x;,x;+1] along the viewing ray, 
which is given by T(x,x;)S|x, — T(x, Xi+1)S)x,,,- By defini- 
tion the points x; are on the boundaries of the shadow volume 
of the terrain. Hence they can be found with a shadow vol- 
ume algorithm such as [HHLH0S5]. This algorithm extrudes 
the silhouette edges of objects, as seen from the light; it also 
projects these objects on the near plane to get correct re- 
sults despite clipping. However these algorithms also gen- 
erate many surfaces that do not correspond to a boundary 
between light and shadow (see Figure 5). These false bound- 
aries must be ignored when computing the inscattered light, 
otherwise a wrong result is obtained. Unfortunately, detect- 
ing them is a non local operation that is not adapted to GPU 
(requiring, e.g., the use of multiple passes, or list structures). 


Our solution is to use the shadow volume algorithm to 
compute the total length / of the shadowed segments, and 
to replace them with a single segment of this length at the 
“ground” end of the ray (see Figure 5). The false boundaries 
still cause problem, i.e., an overestimation of /. Here how- 
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Figure 6: Validation. The sky luminance § in fisheye view for sev- 
eral sun zenith angles, in color, and relatively to the zenith lumi- 
nance. With & = 0.1, Bj, = 21075 m—!, By,/B§, = 0.9, g = 0.76 
and Hy = 1.2 km we get the CIE clear sky model, fitted from actual 
measurements (source [ZWP07]). 


ever / can be clamped to the distance between the nearest 
and farthest faces of the shadow volume. This gives the cor- 
rect result in most cases, and an approximate value in the 
others. Our detailed algorithm is the following. We associate 
with each pixel 4 values An, Az, Zmin, Zmax initialized to 0, 0, 
co, 0. In a first step we decrement (resp. increment) An by 
1 and Az by the fragment depth z, and update Znin and Zmax 
with z, for each front (resp. back) face of the shadow surface. 
In a second step we use (see Figure 5) 


l = clamp(Az — An.Zground' 0, Zground — Zmin) 


L ~ Lo + RILo] + RIL+] + Six — T(x, xs) (17) 


Speed 
when looking at the ground or, when looking at the sky: 
[= clamp(Az, 0, Zmax) 


L ~ Lo +R[Lo] + RJL] +T(x)S,, (18) 


=x+iv 


6. Implementation, results and discussion 


Precomputations We have implemented the precomputa- 
tion algorithm on GPU, with fragment shaders processing 
the numerical integration. This is not mandatory but it allows 
us to quickly change atmospheric parameters, and it saves 
disk space (indeed 5 scattering orders are computed in 5 sec- 
onds on a NVidia 8800 GTS). We store T(r,u) and E(r,u) in 
64 x 256 and 16 x 64 textures. We store S(ur, Uu, Uu,,Uv) = 
[C+,Cmr] in a 32 x 128 x 32 x 8 table, seen as 8 3D tables 
packed in a single 32 x 128 x 256 RGBA texture (using a 
manual linear interpolation for the 4'" coordinate). Thanks 
to our optimized parameterization our 4D table has a better 
precision and uses less space than the 3D table of [SFE07] 
(8 MB for S with 16 bits floats vs 12 MB for their 128° tex- 
ture). 


Rendering The rendering is done in four passes: 


e we draw the terrain in the depth buffer only; 

e we draw the shadow volume of the terrain into a An, Az, 
Zmin, Zmax texture. For this we use the ADD and MAX blend- 
ing functions, disable depth write, and use a geometry 
shader that extrudes the silhouette edges (as seen from the 
sun). This shader also projects on the near plane along —s 
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the back faces (as seen from the sun) that are between this 
plane and the sun [HHLH0S]; 

e we draw the terrain and the other objects with aerial per- 
spective, as well as the sky, using Equations 17 and 18. 
If there are transparent objects such as clouds, aerial per- 
spective must be computed for each object, before blend- 
ing. We use An to compute occlusion in R[Lo], and / com- 
puted as above to get xs (see Section 5); 

e we finally apply a global tone mapping function. 


Results We run several tests with a height field and a re- 
flectance texture from Nasa’s Earth Observatory [SVS*05]. 
Results are depicted in Figures 8 and 9. As shown in Fig- 
ures 6 and 7 our model can reproduce with good accuracy 
the CIE clear sky model, fitted from actual measurements at 
the ground level [DK02]. Since the sky color and aerial per- 
spective are computed with a few texture fetches per pixel 
(< 10), our algorithm is quite fast. For instance, for the right 
view in Figure 8 in 1024 x 768, we get 125 fps without light 
shafts on a NVidia 8800 GTS. This includes 5 ms for the un- 
shaded terrain, 0.4 ms for the first three terms in Equations 
17 and 18 and 2.6 ms for the remaining terms (including 
1 ms to evaluate the non linear parameterization). We get 
25 fps with light shafts (i.e., the first two rendering passes 
cost a lot, about 32 ms). By comparison, we get 50 fps with 
our reimplementation of [O’NO5], using ten samples per ray 
(to get the same quality for single scattering, without shafts). 


Limitations A limitation of our method is that the aerosol 
properties are assumed constant, depending only on altitude, 
whereas in fact they can greatly change depending on the 
atmospheric conditions [Slo02]. Since our precomputations 
are very fast we can change these properties quickly, but they 
remain uniform. 


7. Conclusion 


We have presented the first real-time method to render the 
sky and the aerial perspective from all viewpoints, with mul- 
tiple scattering, terrain shadows and light shafts, and correct 
variation with all view and sun angles. This method is based 
on minimal simplifying assumptions that allow us to get an 
approximate solution of the rendering equation, in which 
most terms can be precomputed. This method can easily be 
extended to more complex physical models, with more con- 
stituents or more wavelengths. 


As future work we would like to model the effect of 
clouds on the ground illuminance and on the aerial per- 
spective, to remove the clear sky hypothesis. Indeed with 
many clouds the interreflections between the ground and the 
clouds should be taken into account [BNLO06]. And their ef- 
fect on aerial perspective should also be considered. To our 
knowledge, this has never been done. 


The source code of our implementation is available at 
http://evasion.inrialpes.fr/~Eric.Bruneton/. 
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Figure 7: Validation. The sky luminance relatively to the zenith luminance for several sun zenith and view zenith angles (and null azimuth be- 
tween view and sun directions). Comparison between our model (with & = 0.1, By =2.21075 m—!, Bi Bie = 0.9, g = 0.73 and Hy = 1.2 km) 
and the CIE sky model 12 (based on actual measurements). We note an overestimation near the horizon (view angles near 90 and -90), which 
is also visible in Figure 6. As shown in [ZWPO7] the Preetham model [PSS99] also suffers from this problem, which probably comes from the 
physical models currently used in CG. 


Figure 8: Results. (a), from top to bottom: /SFEO7], single scattering, multiple scattering and photo. With [SFEO7] the shadow does not 
appear due to the missing V parameter. It is too dark with single sattering only. (b) sunset viewed from space. (c) the view used for our 


performance measurements. 
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Figure 9: Results. Our results (no frames) compared with real photographs found on the Web (red frames). The tone mapping may explain 


the sky hue differences on some images compared with the uncalibrated photographs. 
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Views from space for various altitudes and sun positions. 


Views from the ground showing, from left to right, the Earth shadow, the aerial perspective after sunset, sunset, and light shafts at sunrise. 


Aerial perspective during the day, and mountain shadows for various view and sun angles. 
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